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0^ Abstract 

We investigate the local non-linear dynamics of irrotational dust with vanishing mag- 
\Q netic part of the Weyl tensor, Hab- Once coded in the initial conditions, this dynamical 

^5 restriction is respected by the relativistic evolution equations. Thus, the outcome of the 

latter are exact solutions for special initial conditions with Hab = 0, but with no symme- 
tries: they describe inhomogeneous triaxial dynamics generalizing that of a fluid element 
Q_i in a Tolman-Bondi, Kantowski-Sachs or Szekeres geometry. A subset of these solutions 

^ may be seen as (special) perturbations of Friedmann models, in the sense that there are 

trajectories in phase-space that pass arbitrarily close to the isotropic ones. We find that 
c/3 the final fate of ever-expanding configurations is a spherical void, locally corresponding 

^ to a Milne universe. For collapsing configurations we find a whole family of triaxial at- 

tractors, with vanishing local density parameter i7. These attractors locally correspond 
to Kasner vacuum solutions: there is a single physical configuration collapsing to a degen- 
erate pancake, while the generic configuration collapses to a triaxial spindle singularity. 
These silent universe models may provide a fair representation of the universe on super 
horizon scales. Moreover, one might conjecture that the non-local information carried 
by Hab becomes negligible during the late highly non-linear stages of collapse, so that 
the attractors we find may give all of the relevant expansion or collapse configurations of 
irrotational dust. 

Subject headings: cosmology: theory — large-scale structure of the universe — analytical 
methods — gravitation — relativity 
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1 Introduction 



A good deal of the work of theoretical physicists is spent in constructing models for natural 
phenomena out of their theories, this process perhaps eventually being terminated by an 
experimentalist colleague. As a rule, in this practical process some approximation is taken 
in order to reduce the problem to a tractable one, the hope being that the approximation 
considered is sufficiently reasonable that the results produced are still useful, either because 
they directly give a sufficiently good description of some phenomena or - this goal not being 
reached - because through them we nevertheless gain some clue to the next possible steps we 
can take to set up a better model. 

Roughly speaking, we may divide possible approximations in two broad and not necessarily 
disjoint classes.^ 

The first one we may term the class of exact approximations^ where by the deliberate use 
of this contradictory definition we want to indicate all those exact solutions of the theory 
we are using that can be derived under some special assumptions regarding, for example, the 
matter content and/or the boundary conditions. The second class is that of approximations 
for general data: in this we include all those truly approximate solutions of the equations of 
the theory which can be derived by making some ansatz under which the equations noticeably 
simplify, but still accept generic boundary conditions. 

The search for approximate solutions in these two classes outline two strategies: it seems to 
us that these are equally important and therefore complementary, especially in dealing with non 
linear problems such as those posed by Newtonian and relativistic gravity. In looking for exact 
approximations we often find solutions with unexpected behaviours, revealing aspects of the 
full non linear dynamics, while approximations for general data give us a clue to what can be 
expected under general and reasonable circumstances, but only within the range of validity of 
the underlying assumptions. Examples in the first class are spherical or axisymmetric solutions 
of general relativity with matter, often showing singular behaviour that then deserves further 
study. In the second class we have Newtonian and relativistic cosmological perturbations, 
and the Zel'dovich approximation that reasonably describes pancaking in the first stage of the 
non-linear regime of Newtonian gravity. 

Recently, cosmologists's attention has been drawn on relativistic irrotational dust with van- 
ishing magnetic Weyl tensor, Hab = (Matarrese, Pantano, & Saez 1993, Croudace et al. 1994, 
Bertschinger & Jain 1994). Having in mind to follow cosmological perturbations in the non- 
linear regime in the matter dominated era, the assumption of dust is the simplest description 
one can take for the most interesting case of coUisionless matter. The irrotational assumption 
is well justified if we take a broad view, i.e. if we think of the resulting description as valid for 
not too small scales. Then, the kinematical restriction cOab = is an exact approximation, as it 
is well known that an initially irrotational fiow stays irrotational, in the absence of dissipative 
effects. Moreover, according to the infiationary paradigm for the generation of perturbations, 
no initial vector modes are present at the onset of cosmological structure formation. 

^We do not want to give here a formal classification of possible approximations, we merely introduce these 
two classes for the sake of the following discussion; in particular we do not claim that these two classes are also 
exhaustive. For an extended and somehow related discussion, see (Ellis 1993). 
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The assumption of vanishing magnetic Weyl tensor is the most controversial one. The 
question that was posed and left open by Matarrese, Pantano, & Saez (1993) was if the chosen 
setting Hab = was sufficiently wide to accept generic initial data for irrotational dust. More 
recent works (Matarrese, Pantano, & Saez 1994a, b; Bertschinger & Hamilton 1994, Kofman & 
Pogosyan 1994) have indeed proven that generic purely scalar perturbations of a Friedmann- 
Robertson- Walker (FRW hereafter) universe, although giving Hab = at first order (Goode 
1989; Bruni, Dunsby & Ellis 1992), give rise to a non vanishing magnetic Weyl tensor at second 
order, an effect that we may term tidal induction. 

However, if we want to insert the assumption Hab = (having already assumed cOab = 0) in 
one of the two classes above, then this is an exact approximation. This was implicitly assumed 
by Matarrese, Pantano, & Saez (1993) and was in fact proven by Barnes & Rowlingson (1989) 
in an earlier paper for the more general case of a perfect fiuid. In other words, Einstein 
equations admit exact solutions for an irrotational perfect fiuid with vanishing magnetic Weyl 
tensor, cOab = Hab = 0. Some of these spacetimes are completely general from the point of 
view of the algebraic classification of the Weyl tensor, i.e. they are of Petrov type I, and others 
belong to the degenerate type D or 0. The latter are conformally fiat and all known; type D 
spacetimes are those for which the gravitational field is purely Coulombian (Szekeres 1965), 
e.g. as in the case of Scharzschild or Kerr; type I may be characterized by a superposition 
(because the geodesic deviation equation is linear in the curvature) of purely Coulombian and 
transverse fields (Szekeres 1965). As shown by Barnes & Rowlingson (1989), the dust type D 
spacetimes with Hab = are known explicitly, and due to Szekeres (1975a): in the following 
we will refer to them as Szekeres solutions. 

The assumption Hab = is a dynamical restriction (because is a constraint on the tidal 
field) that has two implications: i) the first comes from one of the usual constraint equations 
(the Hab constraint), which restricts the spatial distribution of the shear CTab (because we also 
have iOab = 0); ii) the second one comes from an additional constraint equation that arises 
when we impose Hab = 0, which restricts the spatial variation of the electric tidal field Eab- 

The restriction Hab = may also be regarded as a sort of gauge choice on the curvature 
(rather than on the metric as it is usually considered in general relativity) that perhaps could 
be dubbed "no induction" or "no radiation", as this, for example, implies fiatness if imposed 
on an exact vacuum solution (Petrov type N) for a gravitational wave (Bruni & van der Elst 
1994); a weaker condition, i.e. the vanishing of Hab at infinity, was shown to be physically 
equivalent to a certain boundary condition in supergravity (Hawking 1982, 1983). 

The interesting feature of the cOab = Hab = perfect fiuid spacetimes is that there is an 
orthonormal tetrad associated with the matter 4-velocity which is the simultaneous eigen- 
frame for the shear CTab of the matter fiow and for the electric Weyl tensor Eab- in all but two 
special cases (Barnes & Rowlingson 1989) the vectors of this tetrad are hypersurface orthog- 
onal, so that a coordinate system exists in which the metric CTab and Eab are all diagonal. 
If, in addition to this, the fiow is taken to be geodesic, = (as it is the case for dust), 
the evolution equations simplify a great deal (Matarrese, Pantano, & Saez 1993), reducing to 
only six ordinary differential equations. Thus, assuming the constraint equations are satisfied 
by the initial data, the following evolution of each fiuid element is no more infiuenced by the 
environment or, in other words, it proceeds as that of a separate universe, which therefore may 
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be dubbed silent universe (Matarrese, Pantano, & Saez 1994a, b). 

This will be exactly the point of view taken here: assuming the constraint equations are 
satisfied, we will investigate the local dynamics of silent models in complete generality, mainly 
using the theory of dynamical systems (e.g. Arrowsmith & Place 1982, Arnol'd 1992) but 
always comparing the outcome of the latter with numerical results. In particular, we will not 
need to assume that the initial conditions are necessarily those arising as linear perturbations 
of a FRW universe, although a subset of such initial conditions is also a subset of the initial 
conditions that are accepted by the equations for silent universes (more precisely, there is 
a subset of perturbed FRW initial conditions that satisfies the constraint equations for the 
p = iOab = Hab = case exactly^ so that these initial conditions are also a subset of all 
the possible initial conditions that satisfy the same constraints). Consistently, we will use a 
set of covariantly defined variables, making no reference to a background FRW model. The 
perturbative point of view can always be recovered in the end, although in general, in doing 
this, there is a problem of gauge choice, basically arising from the fact that in general the 
initial singularity does not occur at the same time in inhomogeneous models, as it happens for 
FRW universes. 

The formalism we use will be briefiy reviewed in Section 2. In Section 3 we will show how 
the introduction of a convenient set of five dimensionless variables and a new time variable 
decouples their evolution equations from that for the expansion (the Raychaudhuri equation). 
After a discussion of some general properties of this novel set of five equations, we will analyze 
in Section 4 the subcase of Szekeres models, arising from the simultaneous degeneracy of 
the shear and electric tidal field eigenvalues. Because of this degeneracy, the Szekeres fiuid 
element may be said locally axisymmetric, although these models admit no Killing vectors. 
These models have been extensively discussed in the literature (e.g. Kramer et al. 1980, and 
reference therein), and retain many of the features of the general case, while only three of 
our variables are needed to discuss their dynamics. Thus, the phase-space for these models 
is 3-dimensional, and can be visualized, revealing general properties that also hold in the 
general 5-dimensional case. In particular, our approach will immediately reveal how there are 
two attracting stationary points for collapsing Szekeres fiuid elements, one corresponding to 
a pancake singularity and the other to a spindle singularity, while the final fate of an ever- 
expanding patch of these spacetimes is to fall into an attracting point representing a spherical 
void (a Milne universe). In Section 5 we will come back to the general triaxial case, and we will 
show that silent models in the collapse phase admit an attracting set which is a closed curve in 
the 5-dimensional phase-space, corresponding to the Kasner (vacuum Bianchi I) solutions of 
general relativity. As it is well known (e.g. Zel'dovich & Novikov 1983, Stephani 1990), there 
is a single pancake in this set, while the generic case is that of the cigar or spindle singularity. 
The pancake is the same as that of the Szekeres models, i.e. is locally axisymmetric, but the 
generic spindle is triaxial. Thus, the final fate of a generic silent fiuid element that stops 
expanding is to collapse to a triaxial spindle singularity, although for perturbed FRW initial 
conditions the final configuration could be almost locally axisymmetric (Bertschinger & Jain 
1994). As it is the case for simple homogeneous models (Zel'dovich & Novikov 1983), the final 
stage of collapse (or the generic initial singularity) is dominated by the curvature, and the 
matter is unimportant, as shown by the vanishing of the local density parameter, 0. Again, 
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the attracting point of ever-expanding configurations is that representing a Milne universe. 

Finally, we will give the asymptotic behaviour of trajectories around the stationary points 
of the phase-space of silent universes. This asymptotic behaviour should be seen as the be- 
haviour of a perturbation around the background solution given by the stationary point itself. 
This analysis will further clarify things like the expansion away from a flat FRW model, the 
expansion toward a spherical void configuration, and the collapse to a pancake or a spindle. 

Conclusions are drawn in Section 6, where we also make some conjectures regarding a 
number of possible cosmological applications of silent universe models. 

Throughout this paper we use units c = SttG = f ; our signature is ( — , +, +, +). 

2 Relativistic dynamics 

In this section we will briefiy summarize a formulation of the relativistic dynamics of collision- 
less matter with zero velocity dispersion (dust) in terms of covariant variables that represent 
observable kinematical and dynamical (curvature) quantities, focusing on the case of irrota- 
tional dust with vanishing magnetic Weyl tensor: p = cOab = Hab = 0. We will also give 
the set of evolution equations for these variables as derived by Matarrese, Pantano, & Saez 
(1993) (see also Barnes & Rowlingson 1989). These equations are a specialization to the case 
p = iOab = Hab = of those presented by Ellis (I97I), where a comprehensive presentation to 
the approach followed here can be found^. 

2.1 Hydrodynamical and gravitational field variables 

Let us consider a perfect fiuid with four- velocity u", normalized to = — I. At each 

spacetime point we may define a projection tensor into the rest space of an observer moving 
with the same four velocity (comoving observer): /i"^ = (/"^ -|- with habu"" = 0. With 
and /i"^ the covariant derivative of any tensorial quantity can be split into a time derivative and 
a spatial derivative. In particular, the spatial part of the derivative Ua-fi of is given by Vab = 
ha'^hb'^Uc-d, with Vabu'^ = 0, while its time part is the acceleration = u^^fcU^, which is a space- 
like vector, ii'^Ua = 0. An overdot denotes convective differentiation with respect to proper time 
t of comoving observers (for a general n-tensor A, Aa^^a2...an = ^ai,a2,...a„;6t^^)- It is standard 
to split the tensor Vab into its trace = its symmetric trace-free part aab = t^(a6) ~ |^a60, 
and its skew symmetric part cOab = V[ab] (the symbol [ stands for skewsymmetrization, ( for 
symmetrization). These are kinematical quantities^ as they determine the relative velocity of 
neighboring fiuid elements. In particular, given a fiuid sphere at an initial time, the expansion 
scalar gives its average volume expansion (0 > 0) or contraction (0 < 0), the shear tensor 
(Tab gives its deformation at fixed volume into an ellipsoid, and the vorticity tensor cOab gives 
its rotation with respect to a locally inertial frame. 

It is possible to give an alternative formulation of general relativity (Ellis I97I) using as 
dynamical variables /?, 0, (Tab^ fj^ab^ and the dynamical quantities (in that they determine the 

^Recently the original 1961 review article by Ehlers has been translated into English by G. F. R. Ellis and 
P. K. S. Dunsby (Ehlers 1993). 
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relative acceleration of fluid elements (Szekeres f 965, Hawking & Ellis f 973) Eab = Cacbd^'^W^ 
and Hab = \T]ac^^C ghbd^'^W^ ■, where Cabcd is the Weyl tensor, i.e. the trace-free part of the 
curvature. The symmetric trace-free tensors Eab and Hab are usually named the electric and 
magnetic part of the Weyl tensor, or simply the electric and magnetic tidal field; they are 
also fiow-orthogonal, EabU^ = ElabU^ = 0. The electric tidal field Eab has a straightforward 
Newtonian analogue (E^^ = cj)^^^ — | (5^^V^^), while there is no counterpart of Hab having an 
independent dynamical role in Newtonian theory (Ellis f97f , Kofman & Pogosyan f994). 

In this formulation, the Einstein equations determine a local algebraic relation between the 
"trace" part of the curvature Rab (the Ricci tensor) and the matter content, as described by the 
energy momentum tensor Tab'- Rab = {Tab — \ QabT). The role of field equations is then played by 
the Ricci identities for u", i.e. Ua-d;c — Ua-c;d = Rabcd'uf' ■, from which the evolution equations for the 
kinematical quantities are derived (as well as a set of constraint equations that must be satisfied 
by the initial data), and by the Bianchi identities in the form C^^'^-^d = i?'''"'^^ — | g''^'^R^\ Thus, 
the Ricci part of the curvature is locally determined by matter through Einstein equations, and 
using these to substitute for Rab in terms of Tab into the Bianchi identities one has a set of 
equations in which the trace-free Weyl part of the curvature is determined non-locally by 
matter. Splitting the Weyl tensor Cabcd into the electric and magnetic parts Eab and Hab 
the Bianchi identities become four equations that can be said Maxwell-like, as they resemble 
Maxwell equations. As usual, the energy and momentum conservation equations T^^^ = 
follow from the contracted Bianchi identities. 

Here we are interested in a perfect fiuid with vanishing pressure, p = 0, and with vanishing 
vorticity, cOab = 0. Two immediate standard results follow from these assumptions: i) since 
p = 0, the acceleration vanishes as a consequence of the momentum conservation equation, 
so that the fiuid fiow is geodesic; ii) since cOab = 0, the fiuid fiow is hypersurface orthogonal, 
i.e. there exist space-like hypersurfaces of which is the normal vector; one can use spatial 
coordinates on these surfaces, and proper time along fiow lines as coordinate time, thus defining 
a comoving synchronous coordinate system. 

However, one can also use a frame rather than a coordinate description, i.e. define com- 
ponents of tensors over a tetrad of vectors. In particular, we can use the orthonormal tetrad 
{M",e|^}, UaC'^ = 0, e'^e^ = Sp (a, [3 = 1,2,3; here greek indices label these vectors) so that 
the tetrad components of tensors are scalars (e.g. Stephani 1990) that are actually measured 
by the comoving observers, as for example the matter density p = Tabu'^u^ ■ 

At this point we can make our third fundamental assumption, i.e. we impose the vanishing 
of the magnetic tidal field Hab = 0. As it was shown by Barnes & Rowlingson (1989) it 
follows from the field equations (the divH equation) that, with this sort of gauge choice on the 
curvature, or dynamical restriction, the eigenframe of the shear (Tab and of the electric tidal 
field Eab coincide. Thus, taking the aligned with the eigenvectors of (Tab and Eab, these two 
quantities may be written as 

3 3 
Eab ^ ^ EfyCfya^ab •) '-^ab ^ ^ '-^a^aa^ab •) (1) 
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with 

where a and E are the magnitudes of the shear and electric tidal field. It was also proven by 
Barnes & Rowlingson (1989) that (in all but two special cases) the spatial tetrad vectors e'^ 
are also hypersurface orthogonal, so that a coordinate basis exists in which the metric is also 
diagonal together with (Tab and Eab- In this case we have 

and the metric may be written as 

3 

ds^ = -dt^ + ^ (4) 

a = l 

where the £„'s give the scaling of lengths in the three directions e'^. Then, a local average scale 
factor £ may be defined, so that the mean expansion rate and the expansion rates in the three 
directions e'^ are given by 

that define the £„'s up to a factor which is constant along each fiow line. Thus, £ is the 
geometric mean of the directional scale factors while |0 is the average of their expansion 
rates: 

^=n4^ -, = \j:f- (6) 

a = l ^ a = l 

With these definitions, we can better see the effect of the shear (Tab- An initially spherical fiuid 
element tends to a fiattened configuration if it has two non-negative shear eigenvalues and a 
strictly negative one, which corresponds to having one direction whose expansion rate is lower 
than the mean local expansion rate |0. It tends to an elongated configuration if it has two 
negative shear eigenvalues and a positive one; in this case two directions have expansion rates 
lower than the local average. When two of the shear eigenvalues are equal, i.e. the shear is 
degenerate, the pancake-like configuration is an oblate spheroid, while the elongated one is a 
prolate spheroid. 

Finally, following a standard terminology (MacCallum 1973, Goode & Wainwright 1982), 
we define singularities as: i) point-like if all three — )• 0; ii) cigar (or spindle) if two of the 
— ^ and the other approaches infinity; Hi) pancake if two of the approach finite numbers 
and the other tends to zero as the singularity is approached. In addition to these, we dub as 
cylinder the special case of two of the — )• and the other approaching a constant value. 

2.2 Silent universes 

Having now introduced all the relevant variables we need to describe the dynamics of irrota- 
tional dust with vanishing magnetic tidal field, p = Uab = E^ab = 0; they are: /?, 0, (Ji, (72, Ei 
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E2. These quantities can be seen as components of a 6-dimensional (6-D) position vector X in 
the phase-space PS6 = {/?, 0, (Ji, (T2, £"1, E2}. It was indeed shown by Matarrese, Pantano, 
& Saez (1993) that the dynamics of irrotational dust with zero magnetic tidal field is described 
by the six first-order ordinary differential equations giving the evolution of these quantities 
along each fiow line. After the time evolution of and the (T„'s has been determined by these 
equations, that of the £„'s will be given by integration of Eq.(5). Thus, the history of each 
fiuid element - being given by ordinary differential equations - proceeds with no infiuence 
from the environment other than those coded in the initial conditions, as in the case of lin- 
ear perturbations of a matter dominated FRW universe (e.g. Ellis & Bruni 1989), or in the 
case of the Zel'dovich approximation in Newtonian theory (Zel'dovich 1970). Since there is no 
communication between neighborhood fiuid fiow lines, we may term these as silent universes 
(Matarrese, Pantano, & Saez 1994a, b). 

We may think of the six equations for our variables as a 6-D fiow of the form X = V{X)^ 
where V at each point X is the tangent to the trajectories in the phase-space PS6. The fiow 
V = V{X) is a non linear function of X only, i.e. V does not depend explicitly on time, so 
that the trajectories do not intersect. Thus, we deal with the following autonomous system: 

p = -Qp, (7) 

= -^0'-2ai2-2aia2-2a2'-^/>, (8) 
2 12 

(7i = - 0-2 (CTI + 0-2) - -CTi^ -- 0(71 - £"1 , (9) 

2 12 

0-2 = - 0-1 (CTI + 0-2) - -0-2^ - - 0(72 - £2 , (10) 
£1 = £1 ((7i-(72)-£2 (fTl + 2(72)-0£l-^/>(7i , (11) 
E2 = £2 (fT2-(7l)-£l ((72 + 2(7i)-0£2-^/O(72 . (12) 

The first of these equations is the matter conservation, the second is Raychaudhuri equation, 
while the other two pairs give the time evolution of the shear and the electric tidal field. Using 
(5) and (6) Eq.(7) gives 

as usual, where M = p*il^ and from now on a subscript * will label quantities at an arbitrary 
time t*. As long as one of the £„'s vanishes, one has a density singularity. In the Newtonian 
theory the case of two of the £„'s contemporarily going to zero is regarded as exceptional, thus, 
the general case is that of Zel'dovich pancaking (e.g. Zel'dovich & Novikov 1983, Shandarin 
& Zel'dovich 1989, Peebles 1993). 

However, the relativistic point of view is different, as what really matters are singularities 
in the curvature. Moreover, for simple homogeneous models the collapse phase is Kasner-like, 
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so that two of the £„'s go to zero together, and spindle-hke singularities are generic, while the 
pancake is exceptional (e.g. Zel'dovich & Novikov 1983, Stephani 1990). 

For the inhomogeneous Szekeres models it is known that a cigar singularity can also occur 
(Goode & Wainwright 1982): however, this is rather obvious: since these models have two 
equal shear eigenvalues, the case for the spindle seems to be as likely as that for the pancake. 
Then, the question is what is the generic collapsing configuration in the triaxial case: here 
we will show that the spindle-like singularity is generic in silent universes, since in the last 
stage of collapse (as well as close to the initial singularity) they tend to have a Kasner-like 
behaviour. 

Although the first work dealing with irrotational fiuids with vanishing magnetic Weyl ten- 
sor is that of Barnes & Rowlingson (1989), these authors never discussed their cosmological 
implications. The first cosmological implementation of these models is due to Matarrese, Pan- 
tano, & Saez (1984), who also showed the existence of spherically symmetric (Tolman-Bondi) 
and planar (Zel'dovich) pancake solutions, these solutions arising in a special case in which (Tab 
and Eab are degenerate, i.e. they have two equal eigenvalues. Later, Croudace et al. (1994) 
discovered an instability of the pancake solution against non-degenerate perturbations, and 
suggested that such an instability could be ascribed to the disregarding of the magnetic tidal 
component. Bertschinger & Jain (1994) soon realized that this instability was actually caused 
by the non-linear tide-shear coupling term in the electric tide evolution equation, which is 
stabilizing for spindle-like collapse but generally destabilizing for pancakes, except for specific 
initial data. Bertschinger & Jain (1994), however, noticed that this result was in contradiction 
with the standard analysis of the Newtonian collapse of isolated ellipsoids (e.g. White & Silk 
1979). Matarrese, Pantano, & Saez (1994a) argued that the non-linear tide-shear coupling 
term is a peculiarity of the relativistic equations, which gives the dominant effect whenever 
Hab is disregarded: by a second-order perturbative calculation (see also Matarrese, Pantano, 
& Saez 1994b, Kojima 1994) they showed that, on scales smaller than the horizon size, where 
the Newtonian approximation should apply, the magnetic tidal tensor cannot be neglected in 
the general case. Only for perturbations on super-horizon scales the Hab = condition would 
apply, leading to the preferential collapse of fiuid elements to spindles, independently of the 
environmental conditions. The problem has been definitely solved by Kofman & Pogosyan 
(1994), who obtained the Newtonian limit of the covariant general relativistic equations by a 
1/c expansion. In particular, they showed that the electric tide evolution equation necessitates 
a calculation at order, in which case a non-zero magnetic tidal field component arises 

as a post-Newtonian effect, related to non-local gravitational interactions. They also reached 
the important conclusion that the magnetic Weyl component implies non-local terms, and in 
the Newtonian limit precisely cancels out the tide-shear coupling term in the tide evolution 
equation. 

Croudace et al. (1994, see also Kasai 1992, 1993, and Salopek, Stewart, & Croudace 1994) 
pointed out that the planar solution was one of those found by Szekeres (1975a), and examined 
the dynamics arising in the limit — )• — oo under the approximation of neglecting linear terms 
in the equations. They also used a set of variables introduced by Matarrese, Pantano, & Saez 
(1993), which are defined as non-linear perturbations with respect to a fiat FRW background. 
This point of view is however rather restrictive, and will not be adopted here: it can always 
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be taken at the end of the analysis. 

As we will show in the next section, it is convenient to introduce a set of new dimensionless 
variables in order to study the dynamics of silent universes. However, some information can 
be already extracted from the system above, before introducing these new variables. 

First, we see that the divergence divV = —5 of V depends just on 0, so that the system 
is dissipative when the fluid element is expanding (0 > 0), and exploding otherwise. From the 
especially simple form Raychaudhuri equation takes for irrotational dust, Fq.(8), we see that 
< 0: thus, once < the collapse of the fluid element is irreversible. The only stationary 
point"^ of the above system is the origin O = {/? = 0, = 0, (Ji = 0, (72 = 0, -Ei = 0, E2 = 0}, 
representing an unstable Minkowski vacuum. Therefore, in the collapse \X\ — )• 00, i.e. the 
magnitude of all the six variables of the above system grows unbounded. 

The system above admits various subsystems describing special subcases. The most obvious 
one is that of the vacuum, p = 0. Another one is that of conformally flat models: these are given 
hy El = E2 = {Hah = is our fundamental assumption), a condition that can be maintained 
only in vacuum (/? = 0), or if the fluid is shear-free (cti = (T2 = 0). Finally, we point out that 
these conformally flat vacuum models, p = Ei = E2 = 0^ correspond to Minkowski spacetime 
in a disguised form. The most interesting subcase is that of the simultaneous degeneracy of 
(Tab and Eab^ and will be discussed in the next section. 

3 Dynamics of silent universes 

In this section we will reformulate the dynamical problem for silent universes in terms of new 
variables that will allow us to simply predict the final fate both of expanding voids and of 
collapsing configurations, as well as the type of the initial singularity. This formulation will 
also make evident a sort of time reversal property of the models. 

3.1 Dimensionless variables 

We start by making the following linear transformation in the phase-space PS6: 

a± = i(ai±a2), ^± = ^ (^1 ± ^2) • (14) 

We then have a new fiow Y = V{Y) in PS6 = {/?, 0, (T_|_, (T_, £"_}, given by the 6-D 
system 

p = -Qp, (15) 

stationary point of a system of differential equations such as (7)-(12), X = V{X), is a point Xs in phase- 
space such that V{Xs) = 0. Such a point represents a special solution of the system, and for our purposes here 
we may say that it can be asymptotically stable, simply stable, unstable, or a saddle. An asymptotically stable 
point will attract generic solutions of the system, i.e. it will be an attractorioi the trajectories in phase-space, 
while trajectories around a simply stable point would go around it, without falling on it (e.g. the potential 
minimum in the mathematical pendulum). A saddle point can be reached if very special initial conditions are 
chosen (a set of measure zero in phase-space), an unstable point usually represents an asymptotically initial 
state, i.e. is a repellerioi the trajectories in phase-space (e.g. Arrowsmith & Place 1982, Arnol'd 1992). 
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(16) 



3 3 + 



(17) 



2 (7+ (7_ - - (7_ - ^_ 



(18) 



E. 



■+ 




(19) 



E. 




(20) 



Again, divV = —5 and the only stationary point is the origin. The interesting subcase is 
now apparent, and given by (T_ = 0, which also imphes = 0: these degenerate (cti = (72 
and El = E2) models are those of Szekeres (1975a). The above system is obviously symmetric 
under the simultaneous change of sign of cr_ and corresponding to the exchange of axes 
1 and 2; it reduces to a 4-D system for cr_ = = 0: the dynamics of these Szekeres models 
will be considered in some detail in the next section. Here, we note that, since this dynamics 
takes place on a 4-D subspace SZ4 of the whole phase-space PS6, trajectories in the latter 
can "go around" SZ4. Even with the above mentioned symmetry, this means that we cannot 
restrict the analysis of the dynamics in PS6 to only a part of it unless we also appropriately 
restrict the initial conditions. Finally, we point out that the degenerate Szekeres models also 
appear under the less obvious restrictions cr_|_ = il^- and E^ = ±^E- in the system above. 
These configurations are replicas of the case cr_ = = which correspond to either CTi = (73 
and El = £"3, or (72 = (73 and E2 = E3. 

In both cases of the two systems (7)-(12) and (15)-(20), the interesting dynamics takes 
place at infinity in PS6. One should therefore introduce a 6-D Poincare sphere and a related 
set of mathematically convenient variables, but this would make the analysis of the dynamics 
rather abstract. Instead, we prefer to introduce a set of new scalar variables 0, E_|_, E_, e_|_, e_ 
directly related to physical observables (0 for example is just the standard density parameter). 
These quantities are dimensionless, and related to the previous variables by 



similar variables have been used in a perturbative context (Goode 1989, Bruni & Piotrkowska 
1994). These relations are obviously singular for = 0, but this is not a fact of major concern: 
it is just a sign of the turn-around = 0, as it is the case for example for in a closed FRW 
model. Thus, these new variables will not be very useful to study the turn-around epoch, but 
they prove to be a good choice to study the collapse and the initial singularities. Moreover, 
the advantage of using them is in the simplifications occurring in the dynamics, which is now 
given by a system of the form 




3 



(7± = E± , £± = e± 0^ ; 
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G = QF{G) , 



(23) 



where G is a position vector in the 5-D phase-space PS5 = {0, E_|_, E_, e_|_, £_}. Now the 
expansion is factored out in the system above, with the only exception of the Raychaudhuri 
equation (22), where we got a 0^ factor. The origin is again a stationary point, but in addition 
to it we now have a whole 5-D stationary hyperplane = 0. However, we have to cut out 
this plane from our analysis, as our new variables 0, E_|_, E_, e_|_, e_ diverge on it. Instead, 
we will now split the analysis of the dynamics of our models in two parts, one for > and 
one for < 0. In order to achieve this splitting, we introduce a new "time" variable 



T = ± Qdt = ±31n. 



(24) 



where we use the plus sign in the above definition when > 0, and the minus sign when < 0. 
This choice ensures that dr / dt > whatever is the sign of 0. Moreover, r is a convenient 
time variable in that r — )• oo in getting close to the collapse or the initial singularity, when 
— )• — oo, as £ — )• 0. Moreover, r — )• oo also signals the complete evacuation of an ever- 
expanding configuration: in such a case — )• 0, but £ — )• oo. Denoting by a prime the 
derivative with respect to r, the evolution equations for our variables for < read 



0' 
0' 

E' 



Q + 6E+2 + 2E_2 + io^ 



^0 (36E+2 - 1 + 12E_2 + 0) , 



E_ 



--E+(l + 6E+)-2E_2--0 

- - 2E+ (3E+ - 1) - 2E_2 - - 
3 + y + J Q 

i-3E+ (4E+-l)-4E_2-io 
i-3E+ (4E+ + l)-4E_2-io 



E_ e_ + -E+ , 

D 



3E_ £+ + -E_ 



6 



(25) 
(26) 
(27) 
(28) 
(29) 
(30) 



We see that, by the introduction of the dimensionless variables 0, Ej- and e± and the time 
T, we have achieved many important simplifications: i) the equations for 0, Ej- and e± do 
not depend on the expansion 0; ii) this means that we have achieved a dimensional reduction^ 
since we can now analyse the dynamics of the 5-D fiow G' = ±F[G) (+ for expansion, — 
for collapse) given by the subsystem (26)-(30); Hi) the change of sign of the expansion 
corresponds now to a "time reversal" r — )• — r, under which the only change in the equations 
above is the sign change of the right hand side; iv) this means that the dynamics is completely 
specular under this time reversal, in the sense that the trajectories in PS5 are exactly the same 
for < and > 0, as only the direction of the tangent to the trajectories is changed for 
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T — )• — T, i.e. — )• —0 <^=> G" — > — G". As we will see, stationary points in PS5 will now appear 
for finite values of 0, E±, and e±: because of i) these stationary points will be the same for 
> and < 0, and because of i) and iv) the Jacobian J = dG'/dG changes sign with 
0: in particular, the eigenvalues of J change sign, so that stable points (attractors) become 
unstable (repellers), and completely unstable points (repellers) become stable (attractors). 
Saddle points of course remain saddles. The Raychaudhuri equation (25) is now merely a 
point dependent clock in PS5 (i.e. the time r (24) depends on the point G in PS5), and it 
can be integrated at stationary points. Because of the factor in the right hand side of 
(25) and considering the definition of r, Eq.(24), 0' will be the same under the time reversal 
T — )• — T <^=> — > —0 at every point G of PS5. This is also true for 0, as it is obvious from 
(22). 

In synthesis, the phase-space during collapse is a mirror image of the phase-space during 
expansion, with the mirror reversing the arrow of the time r. 

3.2 Stationary points 

As we already said, stationary points will now appear for finite (constant) values of , Ej- , e±; 
thus, denoting by a subscript 5 these values, the Raychaudhuri equation (25) at stationary 
points reduces to 

0' = ±/30, /3 = i + 2(3E+| + E_|) + i05, (31) 
3 D 

where here and in what follows the top sign holds for the collapse phase < and the bottom 

sign for the expansion phase > 0. This gives 

= 0, e^^^^-^*^ , (32) 



anc 



e(---') = [1 + /30,(t - t,)]^^ ; (33) 



this relation between r and t clarifies that r — )• 00 when a singularity is approached. For 0, £ 
and the £„'s we have 

® " 1 + /30.(t - t,) ' ^^^^ 

£ = £4i + /30,(t-t,)]^ , (35) 

4 = ia*[i + f3Q*{t-t*)Y" , (36) 

where pa = ^(| + E„) and Ei = E_|_ + E_, E2 = E_|_ — E_ , E3 = — 2E_|_. From these relations we 
see that the expanding or contracting behaviour of each model given by a stationary point is 
monotonic and fixed by the sign of 0*. In dependence of this sign there is a singularity, either 
in the past (initial singularity for 0* > 0) at t — = — (/30*)~^, or in the future (0* < 0) at 
t-t, = (/3|0,|)-i. 
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Before proceeding with the analysis of the more complicated 5-D system (26)-(30), we 
consider the degenerate case (E_ = e_ = 0) of Szekeres models. These will now be described 
by a 3-D flow g' = if{g) (again, + for expansion, — for collapse) in the phase-space SZ3 = 
{0, E_|_, £+}. This flow can be visualized, and this will be instructive also for the most general 
case. 

4 Dynamics of Szekeres models 

That the irrotational dust models with vanishing magnetic tidal field and with degenerate (two 
equal eigenvalues) shear and electric tidal field must be the solutions of Szekeres (1975a) was 
shown by Barnes & Rowlingson (1989). 

Szekeres models generalize Kantowski-Sachs, FRW, and Tolman-Bondi solutions (Szek- 
eres 1975a) and were studied in great detail by many authors because of their interest as 
inhomogeneous solutions of Einstein equations with no symmetries (Kramer et al. 1980; Mac- 
Callum 1993). In particular, Szekeres (1975b) considered the collapse of a subclass with finite 
mass; Bonnor (1976) showed that these spacetimes can be matched to a spherically symmet- 
ric Schwarzschild vacuum solution; Bonnor & Tomimura (1976) studied the time evolution of 
the most cosmologically interesting subclass; Lawitzky (1980) gave their Newtonian analogs; 
Barrow & Silk (1981) considered them in their study of the growth of anisotropic structures in 
the universe; Goode & Wainwright (1982) emphasized how these models can be characterized 
by the growing and decaying modes of dust-filled FRW models, and studied the character of 
their singularities in detail (Goode & Wainwright 1992, and references therein). 

Although Szekeres models are usually divided in two classes, their evolution along fiow 
lines is the same in both of them, as noticed by Goode & Wainwright (1982). Since we are 
precisely interested in this aspect of the models, we will not need to distinguish between the 
two classes. 

4.1 An overview of phase— space 

During the collapsing phase < the dynamics of Szekeres models is given by setting; E_ = 
e_ = in (28), (30). Then, we obtain an autonomous 3-D subsystem in SZ3, g' = —f{g)^ i.e. 








(37) 



E' 



-I- 





e. 



-I- 




(39) 



while, at each point g of SZ3, the Raychaudhuri equation gives r = T[g): 




(40) 
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The divergence divf of the flow / in SZ3 is such that 

divf<0 for > - (l + E+ -42E+2) . (41) 

Thus, for < the flow / is converging everywhere in SZ3, except below the parabola 
given in (41). Remarkably, we see from (41) that divf does not depend on e_|_. From (37) we 
see that the trajectories in SZ3 cannot cross the = plane (because of the factor in the 
right hand side), which corresponds to the fact that only the > section of SZ3 is physically 
meaningful. Also, from (37) we see that 0' < for > 1 — 36 i.e. everywhere (41) is also 
satisfled. Then, from these three facts we may conclude that for < generic trajectories in 
SZ3 converge toward the = plane, i.e. we may expect to flnd stationary points, in particular 
attracting ones, on this plane. 

It is immediate to flnd these points gs with a computer algebra system: they are listed in 
Table 1, giving the coordinate of each point in the phase-space SZ3. With a little more work 
one can easily compute the Jacobian J{gs) = ~[df / dg]{gs) at each of these points, and then 
flnd its eigenvalues (z = 1, 2, 3). These eigenvalues are listed in Table 2, with the values they 
have during collapse, < 0. Points with all < are asymptotically stable, i.e. the generic 
trajectory in SZ3 is attracted by one of these points. Conversely, points with all A^ > are 
unstable. Finally, points with at least two of the A^'s having different sign are saddles, and are 
therefore unstable. 

As previously explained in Section 3, the right hand side of Fqs.(37)-(39) changes sign 
under the time reversal r — )• — r, i.e. in considering the expanding phase, > 0. Then, also 
the Jacobian J{g) changes sign, together with its trace divf and eigenvalues A^. Therefore, 
the inequality (41) is reversed in the expanding phase, and the generic trajectories escape from 
the = plane. Since the A^ have different sign during expansion and collapse, the nature 
of the stationary points also changes. In the last two columns of Table 2 we have therefore 
given the type of each point during collapse and during expansion. To further illustrate the 
behaviour of trajectories in phase-space SZ3, we have plotted in Figure 1 the vector flow / 
above points D II and D V, during the collapse phase, each of these points being at the center 
of the bottom of the plotted boxes. From these flgures it is clearly seen how during collapse 
point D II is repelling, although the flow — / above it is directed downward for 0^ 1; the flow 
— / above point D V is directed toward it, and D V is attracting. 

Figure 2 (a) represents the Poincare disk for the whole = plane during collapse, 
< 0, i.e. the bottom of the boxes in Figure 1, and Figure 2 (b) represents the central part 
of the same plane including all the flve stationary points DII-D VI, for > 0. A comparison 
between the top views in Figure 1 and Figure 2 reveals how, above each point and for small 
enough, the trajectories are basically the same as for = 0. 

4.2 Interpretation of the stationary points 

We now give an interpretation of the stationary points listed in Table 1; again, we remark that 
they are the same for expansion (0 > 0) and collapse (0 < 0). As mentioned in Section 3.2, 

■*This disk is a conformal map of the , £_|_} plane, so that its boundary represents infinity in this plane. 
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in general the stationary points represent models that either expand forever from an initial 
singularity (when we chose > 0), or that collapse from infinity to a future singularity (for 
< 0, i.e. the latter are true time reversal of the former ones). The values of /3 and of the 
exponents p„'s in (36) are given in Table 3. Then, point DI represents a fiat FRW model, 
having zero shear and tidal field. Point D II clearly represents an empty conformally fiat shear- 
free void, locally equivalent to an empty open FRW model (a Milne universe). Point D III is a 
vacuum solution we have not been able to identify in the literature (see below). Using its /3 and 
Pa values, it is easy to check that Point D IV represents the degenerate Kasner model with two 
non-expanding directions, and a pancake singularity. Similarly, point DV represents another 
degenerate Kasner model, with two equally expanding directions and a contracting one, with 
a cigar singularity. Point D VI is the limit of a subclass of Szekeres models (see below, and 
Section 5). Point D VII is clearly unphysical, having < 0. Notice that the points D II, Dili 
and D IV are locally equivalent to Minkowski, (M) in the Tables. In general, configurations 
with E_|_ > (< 0) represent prolate (oblate) spheroidal fiuid elements during collapse (the 
inequalities must be reversed during expansion). 

Now, more interesting than the stationary points themselves is the behaviour of trajectories 
in their neighbourhood in phase-space. For Szekeres models, this is given by the eigenvalues 
of the Jacobian J{gs) = ~[df / dg]{gs) at each of the stationary points gs of the system (37)- 
(39): these eigenvalues are listed in Table 2, for the collapsing phase < 0. Points D I, 
Dili, D VI are saddles, i.e. are unstable in any case, except for special initial conditions. For 
example, during expansion the FRW point, D I, can be approached only by initial conditions 
corresponding to the pure decaying mode of linear perturbation theory (see Section 5.2). Point 
D II is a repeller (an unstable star node) during collapse, and an attractor (a stable star node) 
during expansion: it represents the final fate of ever-expanding voids, locally equivalent to an 
empty open FRW (Milne) universe. Point D III is a saddle in the { E_|_ , e_|_ } plane, and has 
a vanishing eigenvalue for the direction out of this plane: it is a vacuum solution we have 
not been able to identify in the literature, but being a saddle and representing Minkowski, it 
seems of no interest. The Kasner points D IV and DV are repellers (unstable nodes) during 
expansion, and attractors (stable nodes) during collapse. Therefore, point D IV represents the 
final fate of locally axisymmetric pancakes, and point DV that of filaments. Point D VI is 
a saddle, but during expansion it has only one positive eigenvalue: thus, there is a subset of 
trajectories (of measure zero) that fall on it; this point indeed corresponds to the limit t — )• oo 
of Szekeres models P I by Bonnor & Tomimura (1976), i.e. there is a subset of special Szekeres 
models that expand and asymptotically tend to D VI. The unphysical point D VII is unstable. 

Let us look again to Figure 2: the two plots show that the trajectories are the same during 
collapse and expansion, but the directions are reversed: the two attractors D IV and D V in 
the upper figure are unstable in the bottom figure, for > 0, while the origin is now stable; it 
represents the final fate of voids. Note that the e_|_ = axis is a separatrix between the e_|_ > 
and the e_|_ < worlds; since Hah = by assumption, and = on the plane in figure, the 
e+ = line represents Minkowski spacetime: in particular the three points on the line from 
left to right are D IV, D II, and DHL The two points above the e_|_ = axis are D VI (left) 
and D V (right). 

Finally, Figure 3 shows the Poincare disk for the the subcase e_|_ = = 0, i.e. for Minkowski 
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models mentioned above, appearing on the e_ = line in Figure 2. Since only E_|_ varies for 
these models, one can resort to the variables (T_|_ and 0, thus showing the transition from expan- 
sion to collapse. Although the trajectories in the resulting {cr_|_, 0} plane are just Minkowski 
spacetime in disguised form, they give an idea of the general behaviour. 

The central stationary point is Minkowski in its standard form (static and shear free). In 
these conformal maps straight lines stay straight: the vertical (T_|_ = axis is Milne universe, 
point D II (expanding in the upper half, contracting in the lower half); the two lines = 
6(T_|_ and = — 3(T_|_ correspond to points D III and DIV respectively. Trajectories enclosed 
between these two lines in the upper half (lower half) of the figure schematically represent the 
behaviour of ever expanding (contracting) voids. Trajectories outside the two lines represent 
configurations that first expand, then recoUapse. At infinity all the trajectories asymptotically 
approach the two lines. 



4.3 Szekeres trajectories 

In Section 5 we will give the asymptotic behaviour around points D I-D VI of the generic 
trajectories representing triaxial configurations. 

In order to show how the solution of the general system of Eqs.(7)-(12) could be possibly 
afforded, we briefiy report here on the solution of the restricted system for the degenerate 
case (T_ = and = 0, i.e. for Szekeres models. We will basically follow the analysis by 
Goode & Wainwright (1982), but we will only focus on the solution of the time evolution 
equations, without paying attention to the spatial dependence of the various quantities, i.e. 
without dealing with the spatial constraint equations. Note that our notation here will be 
slightly different to that adopted by Goode & Wainwright (1982). 

Up to a scale change which is constant along the fiow lines, we can write 

£^{x,t) = £2{x,t) = S{x,t) , (42) 

and 

£3{x,t) = S{x,t)Z{x,t) , (43) 
which implies p = M/ S^Z^ with M a generally space-dependent integration constant, and 

e.3| + |, ...-if. (44, 

Replacing these expressions into the tide evolution equation, we obtain 

M N 

where N is another (generally space-dependent) integration constant. The evolution equations 
for the expansion scalar and the shear are simplified by introducing the function F = Z — ^ 
which leads to the equation 

F + 2-F = — F, (46) 



S S 
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and 

N 

^ - ^2 • 

The latter equation can be integrated once to give 

2 



s 



2N 



(47) 



(48) 



where k is an integration constant; by properly rescaling S and N one can renormalize k in 
such a way that it can only take the values —1,0, +1. 

It is then clear that the function S is equivalent to the scale factor of matter dominated 
FRW models, in the sense that, for each fluid element, S satisfies the Friedmann equation 
for an open, fiat or closed FRW model, depending on whether k is respectively —1, or +1. 
The function F, on the other hand, has to satisfy Fq.(46), which is the equation for linear 
perturbations around FRW; one has 



F 



(49) 



where /("'"^ and f^~^ represent the growing and decaying modes of linear perturbations around 
FRW and l3^~^\ set the corresponding (generally space-dependent) amplitudes. 

The solutions of Fq.(48) and (46) can be found e.g. in Peebles (1980), we nevertheless 
report their form here for completeness. As far as the scale factor is concerned one has 



S = Nh'{r]), t-to = Nh{r]) 



where 



sinh T] 



1] — sm 1] 



k 
k 
k 



-1 , 
0, 

+1 ; 



(50) 



(51) 



here a prime denotes differentiation with respect to the development angle rj and to is a (gen- 
erally space-dependent) constant. The growing and decaying modes read 



+) 



f [l-(^/2)coth(^/2)] + l 



w'l ' 



6JV 
S 



[l-(7y/2)cot(7y/2)]-l 



k 
k 
k 



-1, 
0, 
+ 1 , 



(52) 



and 



f f^coth(?7/2) 



k 
k 
k 



-1 
0, 
+ 1 



(53) 



Goode & Wainwright (1982) were able to show that, for k = —1 or 0, the final singularity 
in Szekeres models can only be pancake-like and corresponds to Z — )• 0; for k = the 
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final singularity is point -like when /J^"*") = Z?^ correspondence to — > 0), cigar when 

I3(+) < /3(-)/7r (still for S ^ 0) and pancake when /?(+) > /3(-)/7r (for Z ^ 0). 

Finally, let us notice that the planar pancake solution obtained by Matarrese, Pantano, & 
Saez (1993) is recovered in the above formalism, by taking k = = f^~^ = 0. This partic- 
ular trajectory undergoes pancake-like collapse to a shell-crossing singularity, asymptotically 
approaching our attractor D IV. Indeed, one may argue that the only pancake-like collapse 
to D IV (with its replicas) always leads to shell-crossing, contrary to the most generic spindle 
one. 



5 Triaxial dynamics 

In this section we study the stationary points of the triaxial general system (26)-(30), with 
special emphasis on the attracting set, which we will argue to be given by the whole Kasner 
family of vacuum solutions of general relativity. Then, we will give the asymptotic behaviour 
of trajectories in the neighborhood of each of the points DI-DVI, for an arbitrary initial 
condition in this neighborhood in the 5-D phase-space PS5. 

5.1 Kasner attractor set for triaxial collapse 

The triaxial dynamics is similar to that of Szekeres models; again, during collapse trajectories 
are directed toward the = plane, as (26) implies 0' < 0, and 

divG <0 for > 1 - 12 E% (54) 

where = 3E_|_^ + E_^ is the dimensionless shear magnitude. There is however a very 
important difference: in addition to the points D I-D VII (which obviously are stationary 
points also for the general case), and their replicas, there is now a whole family of stationary 
points that form a closed curve in PS5, and which are attractors during collapse. 

The new points are listed in Table 4. Again, these points are the same for expansion and 
collapse. For each listed point T there is a specular one T given by the map E_ — )• — E_, 
e_ — > — e_. More in general there are six distinct replicas (three in T and three in T) for each 
truly triaxial physical configuration, corresponding to the 3! permutations of the principal axes 
of an ellipsoidal fiuid element; each degenerate spheroidal configuration, on the other hand, 
shows up in 3!/2!=3 replicas. 

The eigenvalues of the Jacobian J = —F[Gs) at each stationary point Gs of the system 
(26)-(30) are given in Table 5. The value of /3 in (36), and the exponents are given in 
Table 3. Points TI and T II are replicas of Dili and D VI respectively. Points of the family 
Till are parametrized by —1/3 < E_|_ < 1/3, where the functions in Table 4 read: 

E_(E+) = 1=^1-9J:+\ (55) 
£+(S+) = ^E+(6E+ + l)-i, (56) 



18 



= -^(6S+-l)Vl-9E+^ (57) 

therefore, the family Till is a curve in the 5-D phase-space PS5 (actually, it lies on the 
hyperplane = in PS5) parametrized by E_|_. Let us call A = «4(E_|_) the curve parametrized 
by E_|_, and including both the Till and the Till family. This curve is projected in the 
{E_|_, E_} plane on the ellipse 

3E2 = 9E+2 + 3E_2 = 1 ^ a^ = -Q\ (58) 

3 

This ellipse can obviously be represented in parametric form, using the dimensionless shear 
magnitude E and the angle a in the plane {E_|_ , E_} (see Appendix A). Similarly, the curve 
A in PS5 is projected in the {e_|_ , e_} plane through (56)-(57), and it may be represented in 
parametric form using the dimensionless tidal field magnitude = 3e+^ + and the angle 
/3 in the plane {e+, £_}. Since this curve is parametrized by E_|_ only, /3 is a function of a; 
using [3 = —^a we get for these curves (the subscript yt labels the values on the curve A): 

^+ = 7fS^cos(a), 1 

. with E^ = ^ , (59) 

E_ = Eyt sin(a) , v3 

e+ = -j=eA{/3) cos(/3) , 2 

with £^(;3) = cos(3/3) . (60) 
e_ = e^(/3) sm(/3) , 3V3 

These curves are shown in Figure 4, where that in the {e_|_ , e_} plane is a three-lobe curve 
with period tt (the function eyt(/3) is itself a symmetric three-leaf curve with period |7r). The 
fact that the latter intersect itself is a projection effect: the curve ^ on = in PS5 is closed 
with no self-intersections [A can be visualized in the 3-space { e_|_ , e_ , E_|_ } using E_|_ as 
parameter along the curve). The curve A include all the six possible replicas of each physical 
configuration in the family Till, Till, with degenerate (Szekeres) configurations appearing 
three times (dots in Figure 4, at the intersections with the lines E_ = 0, e_ = 0, E_ = ±3 E_|_ 
and e_ = ±3 e+ (dashed lines in Figure 4. On the curve A^ points T IV, T IV are replicas of 
DVII: all unphysical stationary points are degenerate. Point Till with E_|_ = —1/3 is just 
D IV, and points Till and Till with E_|_ = 1/6 are its replicas: these are the only solutions of 
the family which end up in a pancake singularity. Point Till with E_|_ = 1/3 is just D V, and 
Till, Till, with E_|_ = —1/6, are its replicas; these solutions end up in a cigar singularity. 

During collapse (expansion) points in T III represent fiattened (elongated) ellipsoids for 
— 1/3 < E_|_ < —1/(2^/3) and for < E_|_ < l/(2-\/3), and elongated (fiattened) ellipsoids for 
-1/(2^3) < E+ < and for 1/(2^3) < E+ < 1/3. 

Using in Fqs.(4) and (36) the values of /3 and p„ given in Tab. 3 it is easy to check that 
the whole family T III, T III locally corresponds to the whole set of Kasner vacuum solutions 
of general relativity. As it is well known (Zel'dovich & Novikov 1983, Stephani 1990), there is 
a single pancake singularity in this set (given by D IV and its replicas), while the singularity 
for the generic Kasner model is spindle-like. 
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In Table 5 we give the eigenvalues of the Jacobian J{Gs) at each of the stationary points 
Gs of the general system (26)-(30), for < (again, the A's change sign for > 0). We 
include here also the points D I-D VII of the degenerate case listed in Table 1, in order to show 
their stability properties in the full 5-D phase-space PS5 (cf. Table 2). In particular, the flat 
FRW point D I is again a saddle, and therefore unstable both during expansion and during 
collapse; the Milne point D II is again the attractor for expanding voids. Points TI, T II and 
T IV are equivalent to D III, D VI and D VII respectively, as is made clear by inspection of 
their A's. Also, points D IV and D V are included in the family T III, and their A's are obtained 
by those of Till setting E_|_ = —1/3, E_|_ = 1/6 and E_|_ = —1/6, E_|_ = 1/3 respectively. Note 
that point D II is asymptotically stable in the expanding phase (all A's are negative in this 
case). 

All saddle points are degenerate, and have at least two eigenvalues of different sign; the 
negative ones correspond to eigenvectors in the direction of special trajectories along which 
that point is approached. For instance, in the collapse phase the closed FRW model and a 
special Szekeres model end up in point D I. 

Points of the family Till have A2 = 0, and A^ < (z 7^ 2) during collapse. Thus, the 
outcome of linearization stability analysis is that in the collapsing phase the curve A is an 
attractor. Fach point on A given by a value of E_|_ is asymptotically stable against a perturba- 
tion in a 4-D subspace in PS5, and it is simply stable in the direction of the nearby point given 
by E_|_ + (5E_|_. In other words, the vanishing of one of the eigenvalue indicates an invariance 
in the direction tangent to the attractor curve in phase-space, so that for each point of the 
curve A there should be a 4-D subspace from which this point is reached. This is in fact the 
outcome of any numerical test we have done. Then, this 4-D subspace would be made by 
all those trajectories that end up on that point. Therefore, we conclude that each point on 
the curve A^ locally representing the Kasner set of vacuum homogeneous solutions of general 
relativity, is an attractor for a set of generically triaxial configurations, which is precisely the 
set of trajectories forming the 4-D subspace passing through the point. Thus, the generic 
triaxial configuration that passes through turn-around tends to one point on the curve A^ and 
collapses to a Kasner-like triaxial spindle singularity. Caused by our choice of time variable, 
given by (33), the approach to the singularity is asymptotic in PS5, as expressed by r — )• 00, 
but occurs in a finite proper time t. This approach will be further clarified by the asymp- 
totic analysis of next section. The fact that close to the singularity matter is unimportant, as 
testified by — )• 0, should not surprise. Close to the singularity density and expansion are 
unconnected, so that the density is inhomogeneous, with — )• 00 but — )• 0, while the metric 
can even be homogeneous, as shown for example by Bonnor & Tomimura (1976) for a special 
Szekeres model. 

We have to point out that the discussion above, although proving the attracting nature of 
each single point of the curve A^ is in the following sense incomplete. Indeed, one should be 
able to demonstrate that there is a conserved quantity Q (i.e. Q' = 0), parametrized in the 
same way as the attractor curve A (e.g. by E_|_), that foliates the phase-space PS5. Then, for 
each point Gyt on A there would be a surface passing through it (the 4-D manifold above), and 
the trajectories lying on this surface would asymptotically fall on Gyt (the negative eigenvalues 
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at Gyt corresponding to motion in the surface); the eigenvector of J{Gjs,) corresponding to 
the vanishing eigenvalue, being tangent to the curve, would be transverse (not necessarily 
orthogonal) to these surfaces, and would map the motion from one surface to the next one. 
We stress again that, although we have not been able to show the existence of such a conserved 
quantity numerical tests indicate that it should exist. 



5.2 Asymptotic behaviour around stationary points 

From Eq.(36) one can see that the behaviour of the £„'s is similar to that of Kasner models 
for point D III, while they all expand or contract for D VI (they are trivially the same for D I 
and DII), and all other points are in the Till family, i.e. they are Kasner models. 

We can now give the asymptotic behaviour of the general solutions of the system (26)-(30) 
in the neighborhood of each of these points. This is achieved by using the Jacobian of the 
system at each stationary point, solving the linear problem G' = J{Gs)G with an arbitrary 
initial condition G*, where however G* is supposed to be close to Gs- In doing so, we solve a 
direct problem, i.e. we answer the question of where a trajectory starting in the neighborhood 
of Gs is going to end up to. We leave on a side the inverse problem, i.e. that of determining the 
set of initial data in the neighborhood of Gs from which trajectories would fall on Gs itself. 

For all points except D I, DII and D VI we give here the behaviour for collapse, < 0, 
in terms of the time r; the behaviour in proper time is easily obtained using (33). For DII 
we give the approaching behaviour, i.e. during expansion. Note that the exponents of the 
exp(T — T*) function are the (the eigenvalues of the Jacobian) listed in Tables 2 and 5. 
In what follows 0*, E-|-*,e-|-* are arbitrary initial values for the corresponding variables: the 
closer these constants are chosen to the corresponding values of the variables at the stationary 
point, the closer is the asymptotic solution to a true trajectory falling in that point. 



Point DI /3 = | 0*>O 

= 1 + (0, - l)e3(---') , (61) 

S± = ^(E±,-3e±.)e3(----) + ^(E±, + 2e±,)e-2(----) , (62) 
5 5 

e± = -7(E±,-3e±.)e3(---') + i(E±, + 2e±,)e-2(---') . (63) 
5 5 

Point D I represents a flat FRW universe, and the integration of the linearized system about 
it shows well-known behaviours. It is a 5-D saddle point that repels along the axis, giving 
the known instability of the flat FRW model among curved ones (the flatness problem); then, 
there is a plane containing the two directions Ej- — 3 e± (corresponding to the degeneracy of the 
A = I eigenvalue) from which D I is escaped, and a plane, containing the directions Ej- -\- 2e± 
(for A = — |), from which is approached. The normals to these two planes give the directions 
for the purely decaying and purely growing perturbation modes about a flat matter dominated 
FRW universe. It is in fact straightforward to check, using (33) and /3 = |, that the two modes 
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above go as 1 3 and t ^, as expected. 
Point D II /3 = I 0* > 

= 0,e-3("-"') , (64) 

E± = [E±,-e±,(T-T,)]e-^("-"') , (65) 

e± = e±,e-3("-"') . (66) 

Point D II represents a Milne universe: any volume element that falls in its basin of attrac- 
tion becomes a spherical void. The fact that the fate of ever-expanding configurations is a 
spherical void is a well-known result of Newtonian theory (Icke 1984; see also Bertschinger 
1985), and it is quite interesting that we recover the same result in the fully relativistic context. 

Point D III /3 = I 0* < 

= 0,, (67) 

S+ = -i+ (E+,-e+, + io, + i)e-^(^-^') , (68) 
E_ = (E_, + e_,)e2 -e_,e-2 , (69) 
e+ = (£+* + ^0,)e^(---)-^0, , (70) 

e_ = e_*e-2(-— ) (71) 

This saddle point is a particular Szekeres solution, which belongs to the subset of trajectories 
obtained by setting Af = = in the solutions of Section 4.3; it corresponds to Minkowsky 
space-time in a disguised form. In order not to escape away from it, one has to start with spe- 
cial initial conditions: E_ * = — e_ * and e+ * = — 10*; further, 0* = is required to fall in D III. 

Point D IV /3 = 1 0* < 

= 0,e-("-"') , (72) 

S- = (E_,-£_Je-(---) + £_, , (74) 
£+ = + ^0.) e-2(^--) - ^O.e-(^--) , (75) 

e_ = e_, . (76) 
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This point represents the degenerate (Szekeres) pancake. Starting with generic initial con- 
ditions in its neighborhood, the trajectory ends up in a nearby triaxial point of the family 
Till, unless the special initial condition e_* = is chosen, for which the trajectory falls onto 
D IV. This result, as well as the asymptotic behaviour of e± near the singularity, is in com- 
plete agreement with the numerical analysis by Croudace et al. (1994), who first noticed the 
instability of the planar pancake against general (i.e. non-degenerate) initial perturbations of 
the shear and electric tidal field. Taking e_ * as small expansion parameter, the nearby T III 
point is characterized by E_|_ ~ — | -|- |e^^ and e_|_ ~ l^-*? while E_ ~ e_e_*, thus at first 
order in e_ * the difference between D IV and the nearby T III point is seen only in E_ and e_, 
as shown above. 

Point D V /3 = 1 0* < 

= 0,e-("-"') , (77) 
S+ = i-io,e-(^--)+ (e+,-^E+, + ^0,)e-t(^--) (78) 

E_ = ^(2E_, + £_,)- i(E_, + 3e_.)e-t(---') , (80) 
o 5 



2 

9 6 



lo^e-(---) + ^ - ^ E+, + ^ 0,) g-K^"-) (81) 



+ - e+* - ^ + e-t(---) , (82) 

e- = -i(2E_, + £_,) + ^(E_, + 3e_,)e-t(-—) . (83) 
5 5 

This point represents the degenerate spindle. The behaviour in its neighborhood is very similar 
to that around D IV. The generic trajectory falls in a nearby point Till, at first order seen 
only in E_ and e_, and characterized by E_* = — 3e_*. The special initial condition required 
to fall onto D V is e_* = — 2E_*, explicitly showing that this point is an attractor even for (a 
set of measure zero of) triaxial configurations. 

Point D VI /3 = I 0* > 

= 0,e-4(-— ) , (84) 
^ -l+gE,.4.,.-|0. + l)e-|(^-^^) (85) 

+ S+ . - ^ . - I 0. + ^) ei(--) + I 0.e-i(--) , (86) 
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T — T* 



'16 e_ * — 3 E_ *) sin 
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T — T* 
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0, + — ) e-t(---) 



+ 



5 5 95 



25 
672 



72 



e_ * cos 
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e-^(^-^-^(§7) 

(88) 
(89) 

e-il(^-^')(90) 



This saddle point is also a particular Szekeres solution, which belongs to the subset ol trajec- 
tories obtained by setting Af = in the solutions of Section 4.3. It corresponds to a particular 
model obtained by Bonnor & Tomimura (1976) as the t — )• oo limit of a special set (PI in 
their paper) of ever-expanding Szekeres models. It is indeed seen from (85) and (88) above 
that there are special initial conditions for which the only growing mode, given by A2 = | (in 
Table 5 the signs are given for < 0), is suppressed: for this special set, therefore, D VI is an 
attractor. 



6 Conclusions 

In this work we have carried out a detailed study of the non-linear dynamics of irrotational 
dust with vanishing magnetic tidal field. This type of dynamics is completely described by six 
first-order quasi-linear ordinary differential equations to which the set of partial differential 
equations given by Ellis (1971) reduces in the case p = cOab = Hab = 0. From the point of 
view of the theory of partial differential equations this means that, under the above dynamical 
restrictions, the initial value problem for the general equations is automatically reduced to 
the characteristic initial value problem, where the only surviving characteristics are the fiuid 
fiow-lines. Thus, for these models all the information regarding the environment of each fiuid 
element is that coded in the data on the initial time surface: no sound or gravitational waves 
are allowed to exchange information between neighboring fiuid elements after the initial time, 
so that these models may be termed silent universes. However, coded in the initial value for 
Eab at each fiuid element location, there is always an interaction, even beyond the scale of the 
particle horizon, given by the Coulombian (type D) field generated by distant matter: this 
interaction is ultimately responsible for the spindle-like nature of the singularities we have 
found. The difference with Newtonian cosmology is precisely in this, i.e. in the arbitrariness 
of the boundary conditions for Eab in the Newtonian case (Ellis 1971; Zel'dovich & Novikov 
1983). With appropriate boundary conditions, one can always find a Newtonian solution 
corresponding to the relativistic one, but the vice versa is not necessarily true (Ellis 1971, 
Matarrese, Pantano, & Saez 1994b). 

Although the problem we have studied here is by itself an important and fascinating appli- 
cation of relativistic cosmology, it is tempting to try to interpret our silent universe models as 
an approximation, useful for more general situations, i.e. beyond those particular cases which 
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exactly solve all the constraints required by Hab = 0. In particular, one would be interested 
to know if such a formalism can be applied to the non-linear evolution of scalar perturbations 
in a matter-dominated FRW universe, i.e. to the general problem of gravitational instability 
of coUisionless matter (e.g. Sahni & Coles 1994, for an updated review on the subject). This 
was indeed the question originally posed and left open by Matarrese, Pantano, & Saez (1993). 

The obvious question "why should one resort to general relativity in dealing with structure 
formation in the universe?" , can be answered as follows. It is well-known that the standard 
Newtonian approximation has two main limitations (e.g. Peebles 1980): it cannot be applied 
on scales comparable or larger than the universe horizon size, and it cannot be applied during 
the phases of highly non-linear collapse of a perturbation, when the gravitational interaction 
becomes too strong and/or relativistic motions are produced. These two limitations of Newto- 
nian theory lead us to speculate on two possible applications of our formalism: the description 
of our universe on super-horizon scales, and the study of the highly non-linear collapse of 
cosmological perturbations. 

The first application is rather obvious. On scales larger than the Hubble radius, no causal 
communication is possible, so each patch of the universe on that scale is expected to evolve 
independently of the surrounding ones, gravitational radiation and sound waves are simply 
too slow to carry any useful signal. This can be understood as the "ultra-relativistic" limit of 
Einstein equations. This was in fact proven by Matarrese, Pantano, & Saez (1994a, b) at second 
order in perturbation theory, but is very likely to be valid at any order. What emerges is a 
new picture of the universe on ultra-large scales, with large patches of the universe evolving 
non-linearly according to our equations: only a small set among these will be approximately 
FRW, while most of them will either expand forever to end up in a spherical void, or undergo 
collapse toward some Kasner-like configuration. 

The second, more speculative application of the present formalism is related to the late 
phases of expansion, or collapse after turn-around, of scalar perturbations in a FRW back- 
ground. As far as the fate of ever-expanding configurations is concerned, we think that the 
situation is almost completely settled down: our attractor for the expansion phase (D II) is 
the classical spherical void, already discovered in the Newtonian context (Icke 1984; see also 
Bertschinger 1985, Bertschinger & Jain 1994). That during the late, free expansion of an un- 
derdense region the surrounding matter has no dynamical effect should cause no surprise: in 
such a case the gravitational field does not grow too large, while the overall tendency to make 
the perturbation spherical prevents the occurrence of a relevant magnetic component. In such 
a case the Newtonian approximation already gives the right answer to the problem. 

More problematic is the dynamics of collapse. The present status of our understanding of 
the problem leads to the following picture of the history of a scalar perturbation. At early times 
the magnetic component is negligible (Goode 1989; Bruni, Dunsby, & FUis 1992), as it only 
appears as a second order effect (Matarrese, Pantano, & Saez 1994a, b). Later on, however, 
its presence is generally caused by the spatial gradients of the initial perturbation field. At 
this, mildly non-linear stage, the presence of Hah signals a non-zero fiux of gravitational 
information coming from the environment: it tells to the electric tidal field that the state of 
the surrounding matter has changed compared to the initial conditions; as a consequence the 
electric tide modifies the shear field, which, in turn, acts on the fiuid expansion and density, as 
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prescribed by the Raychaudhuri and mass conservation equations. From this moment onwards 
the evolution of a fluid element is a mixture of two competing effects: that determined by the 
imprint of its "private" initial conditions, and that arising later from the influence of the "rest 
of the world". Which of these two effects will be dominant depends on many variables: the 
form of the particular fluid element, the nature of the initial conditions, such as the type of 
coUisionless matter, the statistics of primordial perturbations, and so on. For instance, a recent 
study within Newtonian theory (Fisenstein & Loeb f994) shows that, at least in the case of 
the standard cold dark matter model, the evolution of shapes of fluid elements is primarily 
induced by the external tide, and not by their initial triaxiality. 

What we believe still to be an open problem is what happens much after a fluid element 
has turned around and its evolution has detached from the overall expansion of the universe. 
Maybe, at this stage, neglecting the influence of the surrounding matter, and thus Hab (in the 
absence of pressure gradients), becomes more feasible. In that case, the dynamics we have 
discussed so far should apply, and the flnal fate of the fluid element should be that of ending 
up into some member of the triaxial attractor set for collapse (Till and Till), i.e. in some 
Kasner-like singularity. 

Notice that we are not claiming that there is a one-to-one mapping between the particular 
flnal conflguration and the speciflc initial condition of the collapsing fluid element, because 
this would disregard the influence of the environment in between. Rather, the picture we have 
in mind is that of fluid elements starting their evolution and ending it according to our local 
set of equations, but evolving in the middle in a much more complicated and non-local way. 
There is an important corollary to our conjecture: since the pancake conflguration (D IV and 
its replicas) only attracts exactly degenerate fluid elements, no fluid element (actually only 
a set of measure zero) would undergo pancake collapse, i.e. no shell-crossing will ever occur 
with non-degenerate initial data: the fate of every collapsing element in a pressureless fluid 
should be a triaxial spindle singularity. 

However, if the latter conjecture should prove incorrect, then we should wonder in which 
sense a non-zero Hab component might modify the above picture. We may try to make 
some guess, this problem being related to the more general one of the robustness (Coley & 
Tavakol 1992) of the p = cOab = Hab— = models. First of all, let us remind that all our 
stationary points are also special exact solutions of the most general system, which is obtained 
by accounting for the magnetic Weyl tensor and for non-zero pressure. Then, the magnetic 
tensor can only act as follows: i) it can introduce new non-trivial stationary conflgurations, ii) 
it can stabilize some of our saddle points or repellers. Hi) it can destabilize some of our collapse 
attractors, and iv) it can modify the basin of attraction for some of our stable points, e.g. it 
can act in such a way that the pancake singularity also attracts a subset of triaxial trajectories. 
All these possibilities, except case i)^ can be studied by considering small perturbations around 
our stationary points, with Hab switched on. In some cases we already know the answer: for 
instance, the standard theory of linear perturbations around FRW tells us that point D I is not 
stabilized by tensor perturbations (i.e. by a linear contribution from the Hab component). The 
remaining cases need to be studied, and this will be the aim of our future investigation. In any 
case, there is one conclusion we can anticipate: if one had to flnd that the basin of attraction of 
our planar pancake (point D IV and its replicas) is enlarged by the presence of a non-zero Hab 
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component it would just mean that the shell-crossing singularity of Newtonian theory (e.g. 
Shandarin & Zel'dovich 1989) only occurs in the non-degenerate relativistic problem because 
of the non-linear interaction with the environment: a non trivial conclusion indeed! 
Let us conclude with an overview of few more open questions. 

First, we have attacked here the problem of the time evolution in silent universes, assuming 
that suitable initial data exist that satisfy the constraint equations for these models (see Barnes 
& Rowlingson 1989). Certainly there exists a subset of perturbed FRW initial conditions that 
satisfy these constraints exactly: it will therefore be interesting to further study them, in order 
to see how much restricted they are. More generally, it will be enlightening to completely solve 
the equations for silent universes, generalizing Szekeres ones, as these solutions will very likely 
reveal further interesting features of these models, as well as some intrinsic symmetry they 
must possess (cf. Kramer et al. 1980, MacCallum 1993). 

We have shown that the final fate of generic triaxial configurations that turn around is to 
recoUapse to one of our spindle singularities of the Kasner attractor. Related to this, there 
are at least two open questions we want to point out, namely: if these singularities are real, 
and, being real, if they are naked. The reality of singularities occurring in the collapse of a 
fiuid of dust is a controversial issue, essentially because of the lack of any pressure gradient 
working against gravity (e.g. Rendall 1992). On the other hand, what really matters from the 
point of view of cosmological structure formation is not whether the curvature singularities 
still occur in the presence of pressure gradients, but if the triaxial collapsing configurations 
we have found give any answer to the question if filamentary or pancake-like structures are 
preferred. Regarding the second point, this clearly goes beyond the scopes of the present work. 
Further, the "time reversal" property of our formulation of the dynamical problem for silent 
universes shows that the Kasner set is also the repeller for generic expanding trajectories, 
so that the typical initial singularity is also spindle-like. Related to this are the issues of 
the isotropy of the observed universe, of its initial state, and of the arrow of time. In silent 
universes, as well as in the Szekeres subcase, there are trajectories that spring off point D I, 
representing the Finstein-de Sitter model, showing that it is conceivable to think of an initial 
isotropic singularity: this will be preferred in a quiescent cosmology^ according to which the 
universe was initially highly regular, this fact in turn being motived by Penrose's idea of a link 
between gravitational entropy and the Weyl tensor, selecting the arrow of time (see Goode & 
Wainwright 1992, and reference therein). 

Of course, an alternative to this point of view is provided by the infiationary scenario. We 
will show in a forthcoming article (Bruni, Matarrese, & Pantano 1994) that patches of a silent 
universe with non-zero cosmological constant that do not recoUapse expand toward a locally 
de Sitter universe. 
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A Parametric representation 



In this paper we have used the variables Ej- and e± to describe the shear of the fluid flow and 
the electric tidal fleld. As it was explained in Section 3 this choice is very convenient in order 
to select the subsystem describing the degenerate case E_ = e_ = (the Szekeres models), 
and also in order to obtain the polynomial character of the flux F[G) in Eq.(23) [the right 
hand side on the system (26)-(30)] which has been used to flnd the stationary points of the 
system. This choice, however, somehow privileges the third axis (or the 1-2 plane orthogonal 
to it). Another possible choice of variable is that used by Bertschinger & Jain (1994), which 
is based on using the magnitude of the shear a = •^J \(t ahf^""^ and of the electric tidal fleld 

E = ^J^EahE"-^ ^ and two angles a and /3. We have used these variables only in the parametric 
presentation of the attractor. Here we give the basic deflnitions, and the relations to the work 
of Bertschinger & Jain (1994). Actually, the quantities used by Bertschinger & Jain (1994) 
differ by ours, in that they deflne the velocity fleld using the conformal time, rather than the 
proper time; this difference is however irrelevant to the present discussion. 

Bertschinger & Jain (1994) represent the shear by = ^o'BjQij, with asj < and 

Q,, = diag [ cos(^) , cos(^) , cos(f ) ] . 

One has Q"^ = |, giving cr^ = |((tJ + a"^ + cr|) = + a"^ + (Ti(T2 = 3(7+^ + (T_^ = |cr|j, i.e. 
(Tbj = — -\/3cr. Similar formulas were given for Eab- In (Bertschinger & Jain 1994) it was 
found that (Tbj = ^ — > — oo at the end of the collapse (with 6 the peculiar expansion scalar), 
corresponding to our exact result E = = | (the background FRW part of is completely 
negligible in the phase of advanced collapse, so = ^ at this stage). 

In this work we have used a = — tt, and /3 = so that the relation asj + 2/3bj = Sir 
goes over [3 = The angle a and [3 are the angles in the planes {E_|_, -^E.} and {e_|_ , -^e.} 
respectively (measured counterclockwise), and are used for the parametric plots in Figure 4. 

Using E^ = 3E_|_^ + E_^ and = 3e+^ + e_^, (e = ^) the relations between E± and E , a 
and e± and e , /3 are : 

E+ = ■^S(a) cos(a) , 
E_ = S(a) sin(a) , 

with E and e functions of the point. 



£+ = -i=£(/3) cos(/3) 
e_ = e(/3)sin(/3) 



(91) 
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Figure captions 



Figure 1: The phase-space flow for Szekeres models during collapse < 0, above points 
D II and DV, placed at the center of the bottom plane = of the boxes; arrows length is 
rescaled in each case (Om is 0, ep is e_|_, and Sp is From top to bottom, on the left: a) 

the flow above point D II up to 2 and including point D I, edge view; b) top view of the 
flow above point D I; c) edge view of the flow above point D II. On the right (top to bottom): 
a) the edge view of the flow above point D V up to 2; b) top view above D V; c) edge view. 

Figure 2: The = (vacuum) plane for Szekeres model. Top: the Poincare disk of the 
plane E_|_, e_|_ (conformal mapping of it) is shown for < 0. Bottom: the central region of the 
plane E_|_, e_|_ including all the flve stationary points is shown for > 0. The trajectories are 
the same, but the directions are reversed: the two stable nodes above are unstable for > 0, 
while the origin is now stable; it represents the flnal fate of voids. Note that the e_|_ = axis 
is a separatrix between the e_|_ > and the e_|_ < worlds; since Hab = by assumption, and 
= on the plane in flgure, the e_|_ = line represents Minkowski spacetime: in particular 
the three points on the line from left to right are D IV, D II, and D III. The two points above 
the e_|_ = axis are D VI (left) and DV (right). 

Figure 3: The Poincare disk of the plane 0, (T_|_ corresponding the the e_|_ = axis in Figure 
2. The behaviour of these Minkowski trajectories simply illustrate that of more complicated 
situations. The central stationary point is Minkowski spacetime in its usual form (static and 
shear free). In these conformal maps straight lines stay straight: the vertical (T_|_ = axis is 
Milne universe, point D II (expanding in the upper half, contracting in the lower half); the two 
straight lines = 6(T_|_ and = — 3(T_|_ correspond to points D III and DIV respectively. 

Figure 4: The whole attractor family {Till, Till} as it shows up in the {E_|_,E_} plane 
(upper left) and in the {e_|_,e_} plane (upper right). The dashed lines are E_/E_|_ = 0, ±3 
and £-/£+ = 0, ±3, corresponding to degenerate conflgurations, and divide the planes in six 
physically equivalent sectors. Dots represent degenerate points of the attractor. The spindle 
is the E_|_ = I dots, with its two replicas at E_|_ = — |; the other three points are the pancake. 
Since this latter is conformally flat, it appears at the origin of the {e_|_ ,e_ } plane, where the 
three replicas coincide; the other three points in this plane are the spindle. Lower flgures show 
e+ = e+(E+) and e_ = e_(E_). 
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Table 1: Stationary points for degenerate (D) models 



Point 







e+ 


Collapse 


Expansion 


Model 


DI 


1 








spherical 


spherical 


Flat FEW 


DII 











spherical 


spherical 


Milne (M) 


Dili 





1/6 





prolate 


oblate 


Szekeres (M) 


DIV 





-1/3 





oblate 


prolate 


Kasner (M) 


DV 





1/3 


2/9 


prolate 


oblate 


Kasner 


DVI 





-1/12 


1/32 


oblate 


prolate 


Szekeres 


DVII 


-3 


-1/3 


1/6 






unphysical 



Table 2: Eigenvalues of the Jacobian J{gs^ < 0) for degenerate (D) models, and point type 
for < and > 



Point 


Ai 


A2 


A3 


Collapse 


Expansion 


DI 


-1/3 


-1/3 


1/2 


saddle 


saddle 


DII 


1/3 


1/3 


1/3 


repeller 


attractor 


Dili 





-1/2 


1/2 


saddle 


saddle 


DIV 


-1 


-1 


-2 


attractor 


repeller 


DV 


-1 


-5/3 


-2/3 


attractor 


repeller 


DVI 


1/4 


-1/4 


5/8 


saddle 


saddle 


DVII 


-1 


-1/2 


1 


saddle 


saddle 
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Table 3: Exponents and singularity type for the stationary solutions 



Point 


/3 


Pi 


P2 


P3 


Singularity 


DI 


1/2 


2/3 


2/3 


2/3 


point-like 


DII 


1/3 


1 


1 


1 


point-like 


Dili 


1/2 


1 


1 





cylinder 


DIV 


1 








1 


pancake 


DV 


1 


2/3 


2/3 


-1/3 


spindle 


DVI 


3/8 


2/3 


2/3 


4/3 


point-like 


Till 


1 


T,Pa 


= TpI = i 


spindle or pancake 



Table 4: Stationary points for the general triaxial (T) system 



Point 







E_ 


e+ 


e_ 


Collapse 


Expansion 


Model 


TI 





-1/12 


1/4 








elongated 


flattened 


Szekeres (M) 


TII 





1/24 


1/8 


-1/64 


-3/64 


flattened 


elongated 


Szekeres 


Till 





s+ 


E_(E+) 




£_(E+) 


flattened or elongated 


Kasner 


TIV 


-3 


1/6 


1/2 


-1/12 


-1/4 






unphysical 
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Table 5: Eigenvalues of the Jacobian J[vecGs^ < 0) for the general triaxial (T) models, and point type 



Point 


Ai 


A2 


A3 


A4 


A5 


Collapse 


Expansion 


r> T 


-1/6 






-1/6 


-i/d 


saddle 


saddle 


DII 


1/3 


1/3 


1/3 


1/3 


1/3 


repeller 


attractor 


Dili 





-1/2 


-1/2 


1/2 


1/2 


saddle 


saddle 


DIV 


-1 





-2 


-1 


-1 


attractor 


repeller 


DV 


-1 





-2/3 


-5/3 


-5/3 


attractor 


repeller 


DVI 


1/4 


-1/4 


5/8 


5/16 - ^^/T5/16 


5/16 + zVT5/16 


saddle 


saddle 


DVII 


-1 


-1/2 


1 


1/4 - z7/4 


1/4 + Z7/4 


saddle 


saddle 


TI 





-1/2 


-1/2 


1/2 


1/2 


saddle 


saddle 


TII 


1/4 


-1/4 


5/8 


5/16 - ^^/T5/16 


5/16 + zVT5/16 


saddle 


saddle 


Till 


-1 





-4/3 + 2E+ 


-4/3 - E+ - -^1 - 9E+ W3 


-4/3 - E+ + ^1 - 9E+VV3 


attractor 


repeller 


TIV 


-1 


-1/2 


1 


1/4 - ^^/^/4: 


l/4 + zV7/4 


saddle 


saddle 
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